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Abstract 

O , 

We consider non-gravitational heating effects on galaxy clusters on the basis of 
the Monte-Carlo modeling of merging trees of dark matter halos combined with the 
thermal evolution of gas inside each halo. Under the assumption of hydrostatic equi- 
qq \ librium and the isothermal gas profiles, our model takes account of the metallicity 



evolution, metallicity-dependent cooling of gas, supernova energy feedback, and heat- 



ing due to jets of radio galaxies in a consistent manner. The observed properties 
of galaxy clusters can be explained in models with higher non-gravitational heating 

Qh 1 efficiency than that in the conventional model. Possibilities include jet heating by 

the Fanaroff-Riley Type II radio galaxies, and the enhanced star formation efficiency 

"55 ! and/or supernova energy feedback, especially at high redshifts. 
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1. Introduction 

Energy feedback plays a vital role in a variety of different phenomena and scales of the 
universe. Supernova explosion produces an overdense shock shell in the surrounding interstellar 
medium which triggers the subsequent star formation. The Gunn-Peterson test in the quasar 
spectra has revealed that our universe was reionized at high redshifts (z ~ 6). The recent 
WMAP result indeed indicated that the reionization epoch may be even earlier than previously 
thought, z— 17 ±5 (Spergel et al. 2003). This implies that the energy feedback from first 
cosmological objects, or non- gravitational heating of the universe, was much stronger than the 
conventional model predictions. 

Non-gravitational heating is also believed to have had a significant influence on the scales 

1 



of galaxy clusters. This is clearly illustrated by the well-known inconsistency of the observed X- 
ray luminosity-temperature (Lx-T) relation; LxocT 3 (e.g., David et al. 1993; Markevitch 1998; 
Arnaud, Evrard 1999) against the simple self-similar prediction Lx oc T 2 (Kaiser 1986). More 
recently Voit, Bryan (2001) suggested that the effect of cooling is important in reproducing the 
observed Lx-T relation (see also Wu, Xue 2002). The subsequent simulations (e.g., Muanwong 
et al. 2002; Kay et al. 2003; Tornatore et al. 2003), however, indicate that while the Lx-T 
relation can be explained by purely cooling effect, the observed hot gas fraction requires a 
fairly significant amount of non-gravitation heating (e.g., figure 3 of Muanwong et al. 2002 and 
figure 8 of Kay et al. 2003). 

Physical origin of the non-gravitational heating of the intracluster medium (ICM) still 
remains to be understood. A plausible candidate responsible for the heating is the energy 
feedback by supernova explosions before and/or during the formation of galaxy clusters (Evrard, 
Henry 1991; Kaiser 1991). Previous authors (Cavaliere et al. 1999; Balogh et al. 1999; Pen 
1999; Kravtsov, Yepes 2000; Loewenstein 2000; Wu et al. 2000; Bower et al. 2001; Brighenti, 
Mathews 2001), however, concluded that the excess heating energy ~ IkeV per gas particle 
is required to account for the observed Lx-T relation. This amount of energy seems larger 
than the conventional model prediction of the supernova combined with the standard star 
formation history and the initial mass function of stars. Another candidate of the heating 
sources frequently discussed is active galactic nuclei (AGNs). Since the radiation from AGNs 
is ineffective in heating the ICM, one has to look for the kinetic energy input, from radio jets 
for instance, to efficiently thermalize the ICM (Wu et al. 2000). Inoue, Sasaki (2001) claimed 
that jets from radio galaxies can provide sufficient energy input to thermalize the ICM and 
explain the observed Lx-T relation on the basis of simple analytic estimates, despite several 
considerable uncertainties of the abundances and intrinsic properties of such radio galaxies. 

In this paper, we explore extensively the consequences of non-gravitational heating pro- 
cesses on the observable properties of galaxy clusters. For that purpose, we follow the merging 
history of dark matter halos in the Monte-Carlo fashion and then trace the thermal evolution 
of baryonic gas inside these halos. This approach was first applied by Wu et al. (2000) in the 
context of the ICM heating, and our current method improves their modeling in several aspects 
and also approaches the problem in a complementary fashion; (i) we follow the merging tree of 
dark halos of mass down to 1.7 x 1O 7 M from z = 30 to z = so as to resolve all the halos that 
may cool, i.e., whose virial temperature exceeds 10 4 K (Shimizu et al. 2002). (ii) We simultane- 
ously consider the supernova energy and the jets of radio galaxies as heating sources in addition 
to various cooling processes. Thus the thermal evolution of baryonic gas, star formation history 
and the metallicity evolution are solved in a consistent manner, (iii) We adopt the cooling rate 
which incorporates the metallicity evolution of hot gas in each halo, (iv) The overall strength 
of heating in our model is controlled by the two dimensionless parameters, e SN and e RG . Their 
values are normalized so that e SN = e RG = 1 for our canonical sets of assumptions. Nevertheless 



2 



we survey a wider range of the parameter space taking account of the fact that the nature of 
those sources are poorly understood, especially at high redshifts. (v) We also consider models 
with the enhanced star formation efficiency at high redshifts (z > 7) so as to look for any pos- 
sible implications on the cluster Lx-T relation on the basis of the recent WMAP suggestion of 
the early reionization in the universe (Spergel et al. 2003). 

The rest of the paper is organized as follows. Section 2 briefly describes our method of 
tracing merger trees of dark matter halos. The basic picture of the non-gravitational heating 
processes, supernova feedback and jets of radio galaxies, is presented in section 3. Our model 
to follow thermal and metallicity evolution of baryonic gas inside dark halos is shown in section 
4. We derive constraints on our model parameters from the observed metallicity - temperature 
relation (§5) and Lx-T relation (§6). Section 7 presents further comparison between our pre- 
dictions and the observed properties of galaxy clusters in X-ray band. We also briefly compare 
our results with that of Wu et al. (2000). We consider the enhanced star formation model at 
high redshift in section 8. Finally section 9 is devoted to summary and conclusions. 

Throughout the paper, we adopt a conventional ACDM model with the following set of 
cosmological parameters (e.g., Spergel et al. 2003); the density parameter = 0.3, the cosmo- 
logical constant £1^ = 0.7, the dimensionless Hubble constant h 70 = ^T /(70kms _1 Mpc -1 ) = 1, 
the baryon density parameter Ob = 0.04/if 2 , and the value of the mass fluctuation amplitude 
at 8 h- 1 Mpc, cr 8 = 0.84, where h = 0.7h 70 . 

2. Modeling Merger Trees of Dark Matter Halos 

Our approach begins with constructing realizations of merger histories of dark matter 
halos. Specifically we use the method (Shimizu et al. 2002) that is based on an improved version 
of the algorithm first proposed by Somerville, Kolatt (1999). 

Consider a halo of mass Mi located at redshift z\. The mass of its progenitors, M 2 , 
at a slightly earlier redshift Z2 = z\ + Az(z\) obeys the mass-weighted conditional probability 
function derived in the extended Press-Schechter theory (Bower 1991; Bond et al. 1991): 

where 5 C ^ ~ 3(127r) 2 / 3 /20-D(2:j) (its useful approximate formula may be found in Kitayama, Suto 
1996) is the critical over-density of the mass density field at a redshift of Zi, D(zj) is the linear 
growth rate, and Si = a 2 (Mi) is a mass variance of the density field top-hat smoothed over the 
mass scale Mj. The corresponding number-weighted conditional probability function for M 2 is 
written as 

dN , 1 . M.\ dP , 1 , 

-(M 2 ,z 2 \M l ,z 1 ) = T ±—(M 2 ,z 2 \M 1 ,z 1 ). (2) 



(1) 



dM 2 K ' 1 ' ' M 2 dM 2 
One needs an algorithm to find all the progenitors M\ {% — 1 ~ N) for M\ which satisfy 
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both the mass conservation Mi = Xijli^ an d equation (2). In reality, this cannot be solved 
without the knowledge of the joint conditional probability function for all M\ that is in fact 
not known. Therefore a variety of empirical prescriptions/tricks have been proposed so far 
(Kauffmann, White 1993; Somerville, Kolatt 1999; Sheth, Lemson 1999, for instance). Our 
algorithm basically attempts to select all the progenitors sequentially as long as M\ is larger 
than the resolution mass at z 2 , M res (z2) and the total mass satisfies 

N 

Y J M\<M X - AM acc (< M res ). (3) 
i=i 

The value of the minimal mass M Tes (z 2 ) is chosen so that its virial temperature is 10 4 K below 
which the gas cooling rate becomes substantially low. This mass resolution is important so as 
not to underestimate the cold gas fraction even at high redshifts. Those halos smaller than 
M res (z 2 ) are not separately treated in the merging tree, but they are collectively accounted for 
as an accretion component. This is why the expression (3) has the additional term defined as 

rM Tes dN 

AM acc (<M res )= / dM 2 M 2 —-(M 2 ,z 2 \M 1 ,z 1 ). (4) 
jo alvl 2 

We need to trace merger trees from z = to sufficiently high redshifts (we choose z = 30) 
when cooling and (non-gravitational) heating of gas are virtually not yet important. This is, 
however, practically impossible because the number of progenitors becomes progressively larger 
at higher redshifts. Therefore we decided to prepare two different sets of merger trees at low 
redshifts (z = to z = 7) and at high redshifts (z = 7 to z = 30), and trace each halo of the 
former at z = 7 by using the tree of its statistical counterpart of the latter. The detail of this 
procedure is described in the Appendix. 

First we construct a merging tree realization starting at z = back to z = 7. Then at 
z = 7, we compute masses of all the existing halos, M hll i oi (z = 7), and begin to trace their 
history down to z = 0. For each halo at the later timestep z{< 7), we search for a halo whose 
mass Mhaio(-2) first exceeds twice the mass of its most massive progenitor at z — 7. If this is 
the case, the halo is assigned its formation epoch as Zf = z. This procedure is repeated for all 
halos at any timestep with respect to the mass of the most massive progenitor at its formation 
epoch. In this way, we assign the sequence of formation epochs for all halos along the merging 
tree realization down to z = 0. 

At its formation epoch, each halo in the merger tree is assumed to obey the density 
profile of dark matter: 

f P(z)5 c (M) 

Phaio(r;M)= (r/r B )(l + r/r„) a vir (5) 

( r > r vir , 

where p(z) = QmPco{^ + z) 3 is the mean density of the universe at z, p c o is the present critical 
density, S C (M) is the characteristic density excess, and r s (M) indicates the scale radius of the 
halo (Navarro et al. 1996). The virial radius r vir is defined according to the spherical collapse 
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model as 

/ W \ 1/3 

and the approximation for the critical overdensity A n i = A n i(f2M,^A) can be found in Kitayama, 
Suto (1996). The two parameters r s and r v ; r are related in terms of the concentration parameter: 

(M \ = r J^MA - J^L ( M \ ~°' 13 m 

C[ " Z) - r s (M,z) l + z\lAxlO u h^M Q ) ' [) 

where the second equality implies its empirical fitting function (Bullock et al. 2001; Oguri et 
al. 2001; Shimizu et al. 2003). Finally we impose the condition that the total mass inside r v ; r 
is equal to M vir which relates 5 C to c as 



3 ln(l + c) -c/(l + c 

We assume that the above profile for each halo does not change until the next formation 
epoch for the corresponding descendant halo. At the next formation epoch, the halo profile 
according to equation (5) is rebuilt with the mass at that epoch. 

3. Physical Models for Non-Gravitational Heating 

In this paper, we focus on two physical models for non-gravitational heating of ICM, 
supernova energy feedback and thermalization of jets of radio galaxies. Their basic pictures are 
briefly described here. 

3. 1 . Energy feedback from supernovae 

To roughly estimate the order of the available energy from supernovae, consider the 
following qualitative picture. If a halo with baryonic gas mass of M gas cools and forms stars of 
mass AM* = /*M gas , the extra heating energy is estimated as r] SN f*M gas E K , where ?7 SN is the 
number of supernova events per unit mass of stars and Ek is a typical supernova kinetic energy. 
The value of r/ SN may be computed from the adopted stellar IMF (initial mass function). In 
the case of the Salpeter IMF dN/d\nM oc M~ x with x = 1.35 over 0.1 M < M < 125 M , one 
obtains 

,125Af rfJV 

/ - — dM 

/ M-—dM 

JoaMq dM 

if only Type II supernovae (corresponding to stars with M >8M & ) are considered as the heating 
source. We adopt E# = 10 51 erg. Then the extra heating energy per gas particle is 
AE SN f) SN AMA 



iVgas M gas /m p 



CsN?7sN/*'Tip-E l K 



5 



Eu 



0.007 Mq\ 



- keV/particle, 



(10) 



10 51 erg J 

where m p is the proton mass. This is approximately the amount of energy that is required to 
account for the Ly^-T relation of clusters as discussed in the previous literature, depending on 
the specific values of /*, ?7 SN , E K , and e SN . 

Admittedly the estimates of these values are fairly uncertain. Recent studies of the IMF 
(Kroupa 2001, 2002) indicate that the x of the IMF is about 1.3 for M > 0.5 M and 0.3 for 
0.08 M Q < M < 0.5 Mq. In this case, ?7 SN = 0.013 M^ 1 , a factor of two larger than the above 
value. Furthermore the IMF is likely to be time- dependent, and ?7 SN (and/or the star formation 
efficiency /*) may be signicantly larger at high redshifts given the WMAP implication for the 
early cosmic reionization. Also the value of E-& may have a large scatter and the efficiency of 
the energy to thermalize the ICM is not clear. Thus our strategy is to adopt the fiducial values 
of those parameters for definiteness, and represent all the uncertainty by the overall amplitude 
parameter e SN . Our fiducial set of model parameters corresponds to e SN = 1, but we explore a 
broad range < e SN < 10 to reflect the uncertainty discussed above. 

3.2. Jets of radio galaxies 

The other important heating source that we consider in this paper is the energy input 
by AGNs. Since the radiation from AGNs are very inefficient in thermalizing the ICM, we 
consider a population of AGNs which have powerful jets, Type II of the Fanaroff-Riley radio 
galaxies (FR II); Type I of the FR galaxies is known to have less powerful jets and we neglect 
the contribution of the latter for simplicity (EnBlin et al. 1997; Valageas, Silk 1999; Enfilin, 
Kaiser 2000; Wu et al. 2000; Inoue, Sasaki 2001; Nath, Roychowdhury 2002). 

To proceed further we need the luminosity function of the FR II radio galaxies. Willott 
et al. (2001) proposed three different evolution models. Specifically we adopt the intermediate 
evolution model (model C) which consists of two terms for low-luminosity and high-luminosity 
radio source populations: 

dn RG (L 151 ,z) dn RG:L dn RG ^ 



dlog 10 Li5i dlog 10 Li 51 + dlog 10 L 15 i 



(11) 



where L151 is luminosity of radio galaxies observed at 151MHz. The two terms are assumed to 
have the following forms: 

-^151 



dn 



RG.L 



d\og w L 



"10 



151 



L 



151 



u 



-"1 



exp 



L, 




(12) 



and 



dn 



RG.H 



dlog 10 L 



151 



'151 



exp 



L 



L 



151 



X < 









exp 






< 












exp 










V z h2 J 



z < z h0 



z > z h0 . 



(13) 
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The original values of the parameters in the above equations are derived for = 0, £7a = 
and # = 50 km s" 1 Mpc" 1 (Willott et al. 2001) which are summarized in table 1. We convert 
their luminosity function to that relevant for our cosmological model (Qm = 0.3, Q\ = 0.7 and 
h 70 = 1) using equation (14) in Willott et al. (2001). 



Table 1. The values of the parameters for luminosity function of radio galaxies. 



n w /h 3 7Q Mpc 3 




W^to 2 WHz" 1 sr- 1 


-in 


k 




8.23 x 1CT 8 


0.586 


1.54 x 10 26 


0.710 


3.48 




riho/Zifo Mpc" 3 


fik 


WVo 2 WHz^sr- 1 


ZhO 


Zhl 




4.80 x 1CT 7 


2.42 


1.25 x 10 27 


2.03 


0.568 


0.956 



It still remains to find the average relation between the radio luminosity L 151 and the 
kinetic power in jet L } of the population of radio galaxies. On the basis of a physical model 
combined with observational constraints of hot spot advance speeds and spectral ages for FR 
II radio galaxies, Willott et al. (1999) derived 

/ r \ 6/7 

L ' = L5Xl045 4 5.1xlO^^WHz-sr- 1 ) ^ S ~' < 14 > 
where f- } is a fudge factor reflecting uncertainties of the physical condition inside the jet 
lobes (e.g., departures from equipartition, proton and low-energy electron content). Willott 
et al. (1999) suggested /j = 1 ~ 20, and we use f- 3 = 10 following other independent studies 
(e.g., Leahy, Gizani 2001; Hardcastle, Worrall 2000; Blundell, Rawlings 2000; see also §3.2 of 
Inoue, Sasaki 2001). Moreover Inoue, Sasaki (2001) showed that equation (14) underpredicts 
the observed jet luminosities (Rawlings 1992) by an order of magnitude (see figure 1 of Inoue, 
Sasaki 2001). 

Therefore our model adopts the relation of 

/ j \ 6/7 

L j = 1.5xl0 47 — , ^"o'ergs- 1 . (15) 

J V 5 - lx 10 27 frfo W Hz" 1 sr" 1 / 70 & v ; 

In order to obtain the total jet energy from Lj, one needs a typical value of the lifetime of 
jets, tiifc. Using a simple model of Falle (1991) combined with the data of 7C Redshift Survey 
and 3CRR sample, Willott et al. (1999) derived a constraint 2 x 10 6 yr < t life < 10 s yr for 
-^151 ^ 5.1 x 10 27 hjQ W Hz" 1 sr -1 . So we adopt t u f e = 10 7 yr for all radio galaxies just for 
simplicity and definiteness. The above set of parameters yields a value of the extra heating 
energy per gas particle for a halo of gas mass M gas hosting one FR II radio galaxy: 

A-Erg _ ^j^life 



7V" gas " RG Mg as /m 



V 



2e *c f 1K !l 7 A { 1^-) { keV/particle. (16) 
\1.5 x 10 47 erg/sy \10 7 jear J \ M gas J 
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Inoue, Sasaki (2001) suggest that the efficiency of the thermalization of the jet kinetic energy 
may be 0.25 ~ 0.4 if the work done by jets stops when the pressures inside and outside of jets 
are equal to each other. Again we introduce a fudge factor, e RG , to represent overall uncertainty 
in our fiducial set of physical parameters. 

Equation (16) indicates that the FR II radio galaxies may be another potential candidate 
for the non-gravitational heating sources of ICM. We emphasize, however, the significant differ- 
ence between the supernova energy feedback and the jet heating; while the former is ubiquitous 
and always accompanies the star formation, the latter is a relatively rare event. The number 
density of the FR II radio galaxies, n RG (z), is ~ 2 x 10~ 7 /?| Mpc _3 at z = [c.f., equation (23)]. 
Thus one FR II per dark matter mass of p(0)/n^{z = 0) ~ 2 x 10 17 h^M e . Over the cosmic 
age of 10 10 yr, we may have 10 10 y r Anfe generations of radio galaxies, and thus we expect a 
dark halo of mass ~2x 10 14 (10 7 yr/ti i f e )/i7 1 M Q hosts one FR II radio galaxy on average. This 
implies that a typical rich cluster has been heated by the jet only for a short duration in the 
past. Thus the jet heating needs to be simulated in a stochastic fashion as in our Monte-Carlo 
approach. In contrast, the amount of the supernova feedback is basically proportional to the 
halo mass, and thus can be taken into account statistically. In other words, the jet heating may 
induce significant variations in the observed Ly^-T relation of galaxy clusters, and even produce 
a population of outliers from the mean relation. 

4. Thermal and metallicity evolution of baryonic gas inside dark matter halos 

Given the above two major scenarios for non-gravitational heating of ICM, we describe 
our detailed models to trace thermal and metallicity evolution of baryonic gas inside dark 
matter halos. 

4-1. Gas cooling 

We assume that the hot gas is isothermal with temperature T gas and in hydrostatic 
equilibrium under the gravitational potential of the dark halo. Since we adopt equation (5) for 
the dark halo profile, the corresponding gas density profile is analytically expressed as follows 
(Suto et al. 1998): 



Pbat(r) = Phot,o exp[-Bf(r/r s )], 



(17) 



where 



B 



2c T vi 



(18) 



ln(l + c)-c/(l + c)Tj 




(19) 



x 



We define the virial temperature T vir as 




(20) 
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where G is the gravitational constant, /im p is a mean molecular weight (/i ~ 0.6), /cb is the 
Boltzmann constant. The amplitude phot,o is computed so as to reproduce the total hot gas 
mass in each halo when integrating equation (17) within r = r v ; r . 

Once the gas profile is specified, we compute the cooling radius r coo i within which the 
gas cools at each redshift z or the cosmic time t(z). The gas cooling time-scale at a radius r is 
given by 



where A(T gas ,Z) is the cooling rate of gas with temperature T gas and metallicity Z, and nn(r) 
is the number density of hydrogen (including both neutral and ionized). We compute the 
evolution of Z for each halo simultaneously, and apply the relevant cooling rate A(T gas ,Z) by 
using tables in Sutherland, Dopita (1993). For physically reasonable profiles, equation (21) is 
a monotonically increasing function of r, and r coo i is computed from the condition: 



While t coo i(z) in the above equation conceptually denotes the elapse time that gas in 
each halo can spend for cooling, its exact definition is fairly uncertain. In the second equality, 
we adopted the definition of Shimizu et al. (2002) based on the following simple picture; when 
a halo of mass Mf forms at the formation redshift Zf, its hot gas is supposed to reach the profile 
[equation (17)] instantaneously. This defines the origin of the cooling time for the halo, and 
T coo i is set to the elapsed cosmic time since Zf. In the subsequent timesteps, we neglect the 
change in the hot gas profile, gas temperature, and metallicity, even if the halo mass M grows 
due to mergers. At each formation epoch of the host halo, its hot gas profile is reset to the 
profile [equation (17)] corresponding to the new mass and temperature, and the origin of r coo i 
is replaced by that at the new formation epoch. 

We apply this procedure for all halos, and the cold gas in each progenitor halo is simply 
accumulated according to the merger trees. For halos that do not have any direct progenitor 
at higher redshifts, we assume that they consist of hot gas only with no excess energy and 
no metals, and set their gas temperature as the virial temperature of the halos. The accreted 
gas [to represent the gas components below our mass resolution, see equation (4)] is always 
assumed to join the hot gas component of the halo without excess energy or metals. 

4-2. Gas heating 

At the formation epoch of each halo (i.e., each mass doubling time to be more specific), 
we compute the increase of the cold gas mass AM co id of the halo since the previous formation 
epochs of all progenitor halos. 

The increased amount of the cold gas mass AM co id is assumed to instantaneously form 
stars of mass AM* = (1 — / r h) AM co id while the rest of it / r hAM co id is reheated by the supernova 
energy and returns to the hot gas of the halo. The returned hot gas carries the metals produced 




3 Phot (r) 
2 fmi p 




(21) 



^cool(r C ool) = T C ool(z) = t(z) ~ t(Z{(M{)). 



(22) 
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by massive stars and pollutes the hot gas associated with the halo. Observationally, metal 
(especially, Fe) abundances of hot gas around NGC 1399, a cD galaxy in the Fornax Cluster 
are estimated as l.lioii ^& within 360 kpc (Ikebe et al. 1992), and 1.5-2.0 Z Q within 20 kpc 
(Buote 2002). Theoretical models of metal ejection from elliptical galaxies by galactic wind 
also suggests a higher than the solar abundance (e.g., David et al. 1990), independently of the 
mass of stars at their final stage. Therefore we assume that the metallicity of the returned hot 
gas, Z cj ect, is equal to 2 Z Q . 

In addition to the supernova heating, we incorporate the energy input by radio galaxies 
as follows. First, we compute the number density, n RG (z), of the radio galaxies at z using 
equations (11), (12), (13), and table 1: 

nno(z)= [ — dn ^ G — (Li 5 i,z) dlog 10 Li5i- (23) 

- / I°gl0 i 151,min «lOg 10 1^151 

We assume that the minimum luminosity of the radio galaxies Li5i, m i n is 1.6 x 
10 25 hjQ W Hz -1 sr" 1 at 151 MHz, which corresponds to the FR I/FR II dichotomy in the 
radio morphology (Bicknell 1995; Laing 1996; Fabian 2001). 

Next, we compute the expected number, iV RG (M halo , z), of the radio galaxies for a halo 
A^haio at z during each time step of the merging tree realization At tTee (z): 

M halo n RG (z) At tree (z) 



N KG (M halo , z) — — — — — — — — : — , (24) 

JM mlB Mn ps [M, z)dM t Me 

where the minimum mass of halos that host radio galaxies is set as M min = 10 12 M . In the above 
equation, n FS (M, z) is the mass function of dark halos for which we use the Press-Schechter 
formula for definiteness. Note that the last factor in equation (24) accounts for the number of 
generations of the radio galaxies during At tree (z). 

Finally, we randomly assign the radio galaxies to each halo in merger trees for z < 7 (we 
ignore the jet heating for z > 7) according to the expected number. As in the above treatment of 
the gas, we sum up the number of radio galaxies at each timestep but input the corresponding 
jet heating only at the formation epoch of the host halo. This procedure gives rise to the 
stochasticity in the heating. The radio luminosity of each radio galaxy is randomly assigned 
according to the observed luminosity function. We denote the amount of total energy input at 
the formation epoch by AE RG . 

4-3. Determining the gas temperature with the non- gravitational heating 

We compute the total amount of excess energy AE gSuS for each halo at its formation 

epoch: 

A-Eg as = -Ecx,prog + AE SN + AE RG - AE cool , (25) 

where -E^prog is the total excess energy already stored in all the progenitor halos, AE SN and 
AE RG indicate the excess energies from supernova feedback and jets of radio galaxies since the 
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last formation epoch of the progenitor halos, and AE coo \ is the loss of the excess energy due to 
radiative cooling. 

The total energy loss AE coo \ is calculated as AE coo \ = Xi ! ee X ,iAM co i d)i , where e eX;i = 
3k-B(T gaSji — T viri j)/(2/im p ) and AM co id,j is the increase of cold gas mass since the previous 
formation epoch of each progenitor halo i. 

The excess energy thermalizes the surrounding ICM as well as exerts work against gravity 
(e.g., Wu et al. 2000). If we neglect the contribution of self-gravity of the gas, we can define 
the total gravitational energy of the g follows: 

-Egrav (-^~gas) / Phot(r)^ halo (r)dV, (26) 

where $haio( r ) is the gravitational potential energy of the halo and given as 

$haio(r) = -AnGr 2 s 5 c p(z) 

for the density profile of equation (5). Then the total energy of gas of temperature T gas is 
written as the sum of the above gravitational energy and the thermal energy thermal (^gas) : 

-Egas(-^gas) -^grav(-^gas) "4" -^thermal (-^gas) • (28) 

Since we assume that the hot gas in each halo is isothermal, ^thermal (^gas) is simply equal to 
3Mhot^BT'gas/(2/imp), where Mhot is the mass of hot gas component in the halo. 

Now we are in a position to determine the gas temperature with the given energy bud- 
get. In reality this is not easy as in the empirical approach that we attempt here, since we 
do not know the efficiency of even gravitational shock heating. Therefore we first adopt a 
conventional assumption in modeling galaxy clusters that the gravitational shock heating is 
sufficiently effective to keep the gas temperature T gas equal to its virial temperature T vir . Then 
we adopt that an additional assumption that the excess energy AE gss is exclusively consumed 
to increase the temperature beyond T vir . Thus we solve the following equation for T gas given 
E gas (T viT ) from gravitational shock heating and the excess energy due to non-gravitational heat- 
ing [equation (25)]: 

-^gas(^gas) = -Egas(^vir) "4" AE gas . (29) 

4-4- Diffuse gas and the metallicity evolution 

With the presence of non-gravitational heating, it occasionally happens that the gas 
temperature becomes much higher than the virial temperature of the halo, and hydrostatic 
equilibrium is unrealistic. In order to take account of such cases properly, we assume that hot 
gas is completely ejected from the host halo if the ratio T gas /T vir exceeds the critical ratio /diffuse- 
The ejected gas is not bound to the halo, but assumed to exist as a diffuse gas component with 
no subsequent thermal and metallicity evolution. The diffuse component joins the bound hot 
gas component again when the descendent of the progenitor halo increases its halo mass so that 



ln(l + r/r s ) 
r/r s 



(27) 
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its virial temperature T^ h exceeds T gas / /diffuse of the diffuse gas. Then the diffuse gas returns 
to the hot gas carrying the total gas energy and the metals acquired at the time when it is 
ejected. While we use the value of /diffuse = 10, the results at z = are insensitive to the choice 
for 10 < /diffuse ^5 100 since the hydrostatic configuration itself [equation (17)] is very shallow 
and effectively approximates a profile close to the diffuse component. 
Finally the mass of metals for each halo is traced according to 

M z = M z>prog + Z cject / rh AM cold - Yl ^AM cokM , (30) 

% 

where Mz )Pr0 g is the total mass of metals in all the progenitors, and Zj denotes the metallicity 
in each progenitor halo i. As explained in §4.2, we set Z e j ect = 2Z@. 

In passing, we summarize in table 2 the parameters used in modeling thermal and 
metallicity evolution of gas. Our model is characterized by the three free parameters / rn , e S N, 
and e RG , which represent the fraction of the reheated gas due to the supernova feedback, and 
the (dimensionless) strengths of the supernova heating and the jet heating that we explore in 
the subsequent sections. 

Table 2. Parameters in the present model. 



symbol 


adopted value 


physical meaning 


^cool 


halo mass doubling time 


time for cooling gas in each host halo 


VSN 


0.007 Mq 1 


the number of supernovae per unit mass of stars formed out of cold gas 


Ekin 


10 51 erg 


kinetic energy of supernova explosion 




2Z 


metallicity of the ejected gas 


h 


equation (15) 


kinetic luminosity of jets of an FR II radio galaxy 


ilife 


10 7 yr 


life time of jets of an FR II radio galaxy 


/diffuse 


10 


a ratio to define the diffuse gas component 


/rh 


(free) 


reheated gas fraction of cold gas to ICM 


1 - /rh 




star formation efficiency 


£sn 


(free) 


efficiency of energy input by supernova feedback 


£rg 


(free) 


efficiency of energy input by jets of radio galaxies 



5. Constraints from the metallicity— temperature relation 

In order to avoid working in the three-dimensional parameter space, we attempt to put 
constraints on / rn from the metallicity-temperature relation of clusters. Since the metallicity 
is produced by the reheated gas from supernova explosion, we can neglect the effect of the 
radio galaxy heating (e RG = 0) except when e SN is very small (< 0.5) as discussed below. We 
first generate many realizations of clusters at z = over a wide range of mass systematically 
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changing the values of / r h and e SN . For each cluster, we compute the gas temperature T gas and 
the metallicity Z according to the procedure described in the last section. Also we compute its 
bolometric luminosity as: 

L bol = 4vr jT* A(T gas ,Z)[n H (r)} 2 r 2 dr. (31) 

From those simulated clusters, we select samples whose gas temperature T gas is higher 
than 2.5 keV and whose bolometric luminosity is larger than 0.1L\, o \(T gas ) at its gas temperature, 
where Lboi(Tgas) is the mean bolometric luminosity of clusters in the observational sample 
compiled by Ikebe et al. (2002). The lower limit on the luminosity is introduced so as to 
incorporate the observational flux limit approximately, but our results are insensitive the value. 
For that selected sample of simulated clusters (iV se i in total), we compute the statistics: 

2 1 V"^ [-^model,i -^obs (^~model,i)] 

Xred - J-j— 2^ 72 ' y 6Z ) 

JV sel i= i a obs,Z 



Zobs(T) = [0.298 - O.OO5(T/keV)]Z , (33) 

where Z^oAci^ and T mo d e i,i are the metallicity and the gas temperature of z-th cluster in the 
sample (i = 1 ~ iV sc i). The mean metallicity-temperature relation Z ^ S (T) is fitted from the ob- 
served data of Fe abundances of Fukazawa et al. (1998). The corresponding standard deviation 
in the observed data is computed as 

1 



Nobs 



(1 



obs,Z iu 

1 »n 



^ [^obs.i — -Zobs^obs^)] 2 , (34) 
"obs i = i 

where iV b s (= 40) is the total number of the sample of clusters in Fukazawa et al. (1998), and 
Z bs,i and T bs,j are the metallicity and the gas temperature of the i-th cluster. We obtain 
fobs.z w 0.063 Z Q . 

Examples of the predicted Z-T gas relation are shown in figure 1. As expected, the 
metallicity is mainly determined by the reheated gas fraction / r h- The contours of Xred on the 
frh ~ £sn plane (for e RG = 0) are plotted in figure 2. Note that e SN represents the efficiency 
of the supernova energy feedback in our model, and the supernova rate ?7 SN is assumed to be 
independent of the value of e SN . This is why our model has metallicity evolution even when 
e SN = as long as the hot gas cools and forms cold gas and stars. The resulting metallicity 
increases as / r h, but decreases as e SN since the star formation rate is suppressed due to the 
stronger energy input from the supernovae. When e SN ^ 3, however, the result becomes almost 
insensitive to e SN and simply determined by the value of / r h alone. In this regime, heating is so 
strong and the suppression of gas cooling is almost saturated. Thus the further increase of the 
energy feedback does not change the result. 

Just for the convenience of the analyses below, we obtained an empirical fit to the / r h - 
e SN relation from figure 2 which reproduces the observational metallicity-temperature relation: 



tan 



yC/rh-0.2) 



(/ rh <0.5). (35) 
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Fig. 1. Mctallicity-temperature relations of clusters at z = 0. The dots with vertical 
and horizontal error-bars indicate the observational data Fukazawa et al. (1998). The 
other symbols with vertical error-bars represent the mean and the corresponding stan- 
dard deviation from 10 realization halos sorted according to their gas temperature. 
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Fig. 2. Constraints on (/ r h, £sn) from the observed mctallicity-temperature relation by Fukazawa et 
al. (1998). The contours of x^ cd = 1, 2 and 3 are plotted. The best-fit empirical relation [equa- 
tion (35)] that reproduces the observational mctallicity-temperature is plotted in dot-dashed line. 

It is encouraging that the resulting range, 0.2 < /rh ^5 0.5, is roughly consistent with other 
models of evolution of elliptical galaxies (David et al. 1990, 1991a, 1991b; Elbaz et al. 1995). 

We also repeat the above analysis for e RG 7^ 0, and do not find significant difference as 
long as e SN > 0.6 because heating by radio galaxies is rather stochastic as noted in §3.2 and 
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does not affect the mean properties of clusters compared with the supernova feedback (see also 
figure 5 below). 

6. Constraints from the luminosity— temperature relation 

Let us move next to deriving constraints on the e SN - e RG plane from the observed Lboi 
- T gas relation of clusters. In doing so, we fix the value of / r h as a function of e SN so as to 
reproduce the metallicity - temperature relation. To be specific, we invert equation (35) and 
adopt 

3 

/rh = f~ arctan(e SN ) + 0.2. (36) 

We would like to confront our simulated cluster samples at z = with the observational 
data combined by Ikebe et al. (2002). The latter is a flux-limited sample of the flux limit Su m = 
2 x 10 -11 erg s _1 cm -2 in the 0.1-2.4 keV band with spectroscopically measured temperature for 
T gas > 2.5 keV. The number of these clusters is iV bs = 52 in total. Since our simulated clusters 
are not assigned the distance, we cannot construct the corresponding flux-limited sample in 
reality. Thus we select those simulated clusters whose gas temperature T gas is higher than 2.5 
keV and whose bolometric luminosity is larger than 0.1Lb o i(T gas ) as in the previous section. 



CD 



10 3 



10 2 



10 



0.1 



.f rh : fixed from Z-T relation 

— e SN =0.0,e RG =0.0 
e SN =0.0,e RG =2.0 

^ = 2.0,^=1.0 

e SN =4.0,e RG =0.0 
e SN =4.0, eRG =2.0 _ tr-- 



-x- 
-O- 




T aaa [keV] 



Fig. 3. Luminosity-temperature relations of clusters at z = 0. The dots with vertical and 
horizontal error-bars indicate observational data from Ikebe et al. (2002). The other sym- 
bols with vertical error-bars represent the mean and the corresponding standard deviation 
from 10 realization halos sorted according to their gas temperature. The dashed line in- 
dicates the 'self-similar' relation, in which all baryons in the halo have T gas = T v - lr . 



For such selected samples of simulated clusters, we compute 
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Xred 



g [log 10 (L modeM /10 4 4 erg s" 1 ] 



1 obs,logL ( -^modelji 



sel j =1 

1.17 + 2.81 log 10 (T/keV), 



(7 



logL 



(37) 
(38) 



^obs.logL (T") = 

where L mo deH and T mo d e H are the bolometric luminosity and the gas temperature of the i-th 



simulated cluster (i = 1 ~ iV" se i). The function -Fobs.iogL^) is our best-fit 'log 10 L bo i - log 10 T gas ' 
relation to the data of Ikebe et al. (2002). The corresponding standard deviation in the observed 
data is computed as 



A ! ,, 



'logL 



E 



log 



J obsA 



10 



10 44 erg s" 1 



^obsJogL (T hs,i 



(39) 



-^obs i=l 

where L ^ S)i and T ^ si are the bolometric luminosity and the gas temperature of the i-th observed 
cluster (i = 1 ~ iV b S )- We obtain cri og L = 0.24. 

In order to perform the fit over a wide dynamic range of the bolometric luminosities, 
we decided to do so in logarithmic scales instead of in linear scales. Thus the variable x^ ed 
does not obey the standard ^-distribution. We use the value to examine simply the degree 
of the goodness-of-fit between our model and the observation, and do not intend to assign any 
statistical significance in a strict sense. 



Constraint from L-T relation 
L: fixed from Z-T relation- 




2 4 6 



Fig. 4. Constraints on (e SN , e RG ) from the observed luminosity-temperature relation by Ikebe et al. (2002). 
The contours of Xrod = b 2, and 3 are plotted. We fix the value of / r h according to equation (36). The 
six filled circles indicate the parameter sets whose statistical properties are examined in detail later (§7). 



Figure 3 shows examples of our model L^-T relation against the observational data. This 
clearly indicates that the Lx-T relation is very sensitive to the efficiency of the non-gravitational 
heating processes as expected. The case of purely gravitational heating (e SN = e RG = 0; solid 
line with solid squares) reproduces the slope of the analytical self-similar prediction (dashed 

16 



line). Our model incorporating gas cooling, however, results in a significant fraction of cold gas 
in all clusters (see figure 8 below), and the amplitude is lower than the latter in this case. The 
change of the Lx-T relation is induced through those of the cold gas fraction and the hot gas 
density profile due to non-gravitational heating. 

Figure 4 plots the contours of Xred 011 6 sn ~~ £ rg plane. As is well known, the supernova 
feedback (or more strictly, any uniformly heating source of the ICM) is needed to account for 
the observed Lx-T relation. The jet of radio galaxies alone cannot provide the required degree 
of heating for clusters. Keeping the constraint in mind, we will examine various statistical 
properties of clusters in the next section. In particular we focus on the six sets of parameters, 
(e SN ,e RG ) = (0.0,0.0), (0.0,1.0), (1.0,0.0), (1.0,1.0), (3.0,0.0), and (2.0,0.5), plotted in figure 4. 
Incidentally Kay et al. (2003) conducted similar studies using cosmological hydrodynamic sim- 
ulations. While we compute the SN gas temperature/entropy from the amount of the heating 
energy [equation (29)], they explicitly assume its temperature. Thus direct comparison with 
their simulations needs caution, but still our reading of their table 1 is that the significant 
amount of SN heating (i.e., their e > 1) is required to simultaneously account for the Lx-T 
relation and the hot gas fraction of ICM. This is consistent with our above successful parameter 
set (e SN ,e RG ) = (3.0,0.0) as we show below. 

7. Statistical properties of simulated clusters 

The metallicity-temperature relations of our simulated clusters in the six models are 
plotted in figure 5. We adopt equation (36) for f rh even in the case of e RG 7^ (right panels). 
Since the values of / r h are chosen so as to reproduce the Z - T gas relation for e RG = 0, the 
agreement in left panels is just by construction. Right panels in figure 5 confirm that the non- 
vanishing e RG does not affect the mean relation as mentioned in §5, but mainly add scatters 
around the mean due to the stochastic nature of heating by radio galaxies in our model. Note, 
however, that it is probably premature to compare the dispersions around the mean even if 
tempting. It is not clear to what extent the observational data reflect the intrinsic scatter 
rather than merely the observational errors. Furthermore our modeling is very simplified and 
ignores various possible sources for the dispersions; we assumed constant ?7 SN , f rh , e SN , e RG , and 
so on which are likely to be dependent on the environment. In this sense, we expect that our 
model predictions systematically underestimate the real scatters. This remark should apply 
also to the other statistical properties discussed in this section. 

Figure 6 shows the degree of the goodness-of-fit of the luminosity-temperature relation 
in detail for the six sets of model parameters. Top panels indicate that the supernova feedback 
is indeed essential in explaining the L\, \-T relation. Middle and bottom panels suggest that the 
jet heating is relatively insignificant if efficiency of the supernova feedback is high. Although 
the Xr ed -statistics in figure 4 may indicate a marginally acceptable fit in the case of (csn^rg) = 
(1.0,0.0), further heating (either by the jet or by much higher star formation efficiency as 
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Fig. 5. Scatter plots of the mctallicity-temperature relation of clusters at z = 0. Our 
model predictions are plotted in crosses, and the observed Fe abundances (Fukazawa et 
al. 1998) are plotted in dots with error-bars. The solid lines indicate the mean ob- 
servational mctallicity-temperature relation and the corresponding ± lcr dispersions. 

empirically proposed in §8 below) is indeed preferred. It is interesting to note that the jet 
heating sometimes produces an X-ray dark cluster. This corresponds to a cluster in the middle- 
right panel of figure 6; the mass and the temperature of the cluster is 7.9 x 1O 13 M , and 2.7 
keV while it is very X-ray faint (Lboi ~ 7 x 10 42 erg s _1 ). This is an illustrative example which 
reflects the stochastic nature of the jet heating. 

Figure 7 shows the results for the mass - temperature relation of the clusters, which are 
in good agreement with the conclusions in figure 6; we need strong heating sources to match 
the observed M - T relation (Finoguenov et al. 2001). While this may be due to the fact that 
the M - T relation is not entirely independent of the L^ Q \-T relation (since our current model 
adopts the specific density and temperature profiles as a function of mass alone), it will provide 
additional insight at least (see also Shimizu et al. 2003). 

The successful heating models also reproduce well the observed hot gas fraction of clus- 
ters. In figure 8, we plot the hot gas mass fraction in halos at z = against their gas temper- 
ature, f hot (T gas ) = M hot /(tt B M hSk \ (T gas )/Q M ). Again the average value of f ho t{T gaB ) is mainly 
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Fig. 6. Scatter plots of the luminosity-temperature relation of clusters at 2 = 0. Our model predic- 
tions are plotted in crosses, and the observational data (Ikebe et al. 2002) are plotted in dots with 
error-bars. The solid lines indicate the mean observational luminosity-temperature relation and the cor- 
responding ± lcr dispersions. The dotted lines represent our selection criteria; T gas > 2.5 keV and 0.1 
times the mean Lboi(T gas ). The dashed line indicates the 'self-similar' relation corresponding to figure 3. 

controlled by the amount of supernova energy feedback, and is fairly insensitive to the value of 
e RG . Note that cooling of ICM is so efficient and without non-gravitational heating most of the 
ICM remains cold. This is another reason why a significant amount of non-gravitational heating 
is required in galaxy clusters (again not entirely independent of the luminosity - temperature 
relation, though). 

Ponman et al. (1999) find that the 'core entropy', S c , of low temperature clusters (T gas ~ 
IkeV) is much larger than the self-similar model prediction, where S c = kBT gas /[n e (r c )] 2 / 3 is 
defined at the core radius r c . This departure from the self-similarity of ICM also supports the 
importance of non-gravitational heating/ cooling processes in the evolution of galaxy clusters; 
when some fraction of gas in the ICM cools, the hot gas density responsible for the entropy is 
reduced and the gas entropy increases. Of course the gas entropy is increased by the temperature 
rise and the resulting spreading of the core region due to non-gravitational heating. In reality, 
both processes are intimately coupled with each other and one has to trace the thermal evolution 
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Fig. 7. Scatter plots of the mass-temperature relation of clusters at z = 0. Our 
model predictions arc plotted in crosses. The dots with vertical and horizontal er- 
ror-bars indicate the observational data points of Finoguenov et al. (2001) after ap- 
propriate correction for the virial mass (c.f., Appendix A in Shimizu ct al. 2003). 

of ICM taking account of heating and cooling in a self-consistent fashion. Although this seems 
a straightforward project for hydrodynamic numerical simulations, it is not easy to reliably 
trace the thermal evolution over the required dynamic range between < 1O 7 M (that can cool 
at z ~ 30 corresponding to the virial temperature of ~ 10 4 K) and ^ 1O 16 M (that encloses 
the region surrounding a rich cluster) given the limited mass resolution of current simulations. 
This is exactly why our current approach plays a complementary role in studying the thermal 
evolution of ICM. 

Figure 9 shows the core 'entropy' against the gas temperature of halos at z = 0. In 
the plot, we adopt a simple definition for the core radius, r c = 0.1r vir , which is approximately 
consistent with most observed clusters. We estimate n e (r) from equation (17) assuming that 
the hot gas is in collisional ionization equilibrium. Despite such crude definitions adopted 
tentatively, our heating models fairly reproduce the observational best fit relation by Ponman 
et al. (2003) and the value of the entropy 'floor' by Lloyd-Davies et al. (2000). 

Finally the excess energy per gas particle, AE gas , is plotted in figure 10 against the 
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Fig. 8. Scatter plots of the hot gas fraction of clusters z = 0. Our model predictions are plotted in 
crosses. The observational data points are taken from Mohr et al. (1999), but are multiplied by Om/^b- 

halo mass at z — 0. As expected, all the successful models (middle and bottom panels) have 
AE gas ~ IkeV/particle. As we have seen in other plots, the jet from radio galaxies increases 
the dispersion. Note that the excess energy is defined as per one hot gas particle; the hot gas 
mass fraction in the top-right panel is merely 30 percent and this is why AE gss ranges up to 
~ 2keV/particle although the total heating energy is less than the other models (middle and 
bottom panels). 

As we mentioned in Introduction, our current approach corresponds to the improved 
modeling of the earlier work by Wu et al. (2000). Thus it is relevant to compare the results 
here. They concluded that the excess energy of ~ (l-3)keV/particle is required to reproduce the 
observed luminosity - temperature relation. This is in good agreement with others as well as our 
results (figure 10). They also argued that this amount of excess energy is unlikely to be provided 
by supernova feedback except in a highly contrived model of galaxy formation, and suggested 
the AGNs as a potentially important heating source. In this respect, the normalization of our 
model parameter, e SN , may be a bit optimistic; they assume the average kinetic energy from 
supernova of 4 x 10 50 e Wu SN erg with their efficiency parameter e WUjSN ~ 0.3, while we adopted 
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Fig. 9. Scatter plots of the core 'entropy' of clusters at z = 0. Our model predictions 

are plotted in crosses. The solid line indicates the entropy floor of 124 h^^ 3 keV cm 2 

from Lloyd-Davies et al. (2000). The dotted line is the observational best-fit relation by 

Ponman et al. (2003). The dashed line represents the self-similar model prediction. 

10 51 e SN erg. Since they did not show plots of the statistical properties of their resulting clusters 
as we did in the present section, further quantitative comparison is difficult, but the above 
overall qualitative conclusions seem to be in reasonable agreement with each other. 

8. Model of higher star formation activities at z > 7 

One reasonable interpretation of the recent WMAP result on the early reionization 
at z = 17 ±5 (Spergel et al. 2003) is that the star formation rate was much higher and/or 
the energy injection from supernovae was stronger; it has been argued that the first objects 
at high redshifts may be preferentially very massive (~ 100 M Q ), and that their birth rate 
may be higher than those in the range 1-2 M (e.g., Bromm et al. 1999; Abel et al. 2000; 
Nakamura, Umemura 2001). Such first objects may turn into hypernovae with the explosion 
energy > 10 52 erg. In order to mimic such plausible scenarios, we consider a model with the 
bimodal supernova feedback efficiency: 
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clus- 



"SN 



(40) 



e SN ,L < 7) 
x e S N, H (z > 7) 

For definiteness we fix that e SN L = 0.3, corresponding to the 30 percent conversion efficiency of 
the supernova energy to the ICM thermalization. Since the efficiency at high z is much more 
uncertain, we explore the four cases; e SN H = 1.0, 2.0, 5.0, and 10.0. We neglect the jet heating, 
i.e., e RG = 0, and adopt that / r h = 0.26 [equation (36) with e SN = 0.3] at all redshifts. 

Fig lire XI shows tlicit tliG luKiiiiositjy - tjGiiipGrcitju.rG rdcitioii of the model clusters for 6gN n rs - / 
5.0 (lower left panel) reproduces the observation reasonably well. This may correspond to an 
original idea of preheating by Evrard, Henry (1991) and Kaiser (1991), and thus the result is 
fairly insensitive to non-gravitational heating at lower redshifts. 

The other relations are shown in figure 12 for e SNjL = 0.3 and e SNH = 5.0. Since we 
did not attempt to tune the parameter sets, we would say that the agreement is satisfactory; 
the inclusion of the jet heating by radio galaxies at low redshifts would certainly improve the 
agreement. These results are impressive because we do not have to consider exceedingly strong 
supernova feedback at low redshifts if the star formation activity at high redshifts is sufficiently 
enhanced. 
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Fig. 11. Same as figure 6 but for models with enhanced star formation activ- 
ity for z > 7. We assume that e RG = and / r h = 0.26 at all rcdshifts. 

9. Summary and conclusions 

We have developed a semi-analytical approach to trace the thermal history of galaxy 
clusters based on the Monte-Carlo modeling of merging trees of dark matter halos. Under the 
assumption of hydrostatic equilibrium and the isothermal gas profiles, we have incorporated the 
metallicity evolution, the metallicity-dependent cooling of gas, the supernova energy feedback, 
and heating due to the jet of radio galaxies in a consistent manner. The latter non-gravitational 
heating processes were characterized by two dimensionless parameters, e SN and e RG , and we 
explored several statistical properties of galaxy clusters over a wide range of the parameter 
space. 

As has been known for a while, we confirmed that a fiducial model of supernova feedback 
alone, i.e., e SN < 1, does not reproduce the observed luminosity - temperature relation of clusters. 
A reasonable agreement can be achieved by enhancing non-gravitational heating in two different 
ways; i) considering additional heating due to the jet of some class of AGNs, notably Type II 
of the Fanaroff- Riley radio galaxies, and ii) adopting somewhat higher star formation efficiency 
and/or supernova energy feedback. The former possibility was first examined seriously by Inoue, 
Sasaki (2001), and the current study basically confirmed their conclusion using a significantly 
improved methodology to trace the thermal history of ICM. The latter idea is particularly 
interesting in the light of the recent WMAP finding of the earlier reionization epoch of the 
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Fig. 12. Scatter plots of the metallicity-temperature relation, the mass-temperature relation, the 
hot gas fraction, and the core 'entropy' for models with enhanced star formation activity 
for z > 7. Our model predictions are plotted in crosses. We assume that e SN L = 0.3 
and e SNiH = 5.0, while the other parameters have the same values as in figure 11. 

universe than previously thought. By increasing the feedback efficiency at high redshifts (e SN ~ 5 
at z > 7, for instance), most of the model predictions of the simulated clusters can be brought 
into agreement with the observational data. 

So far we discussed several properties of clusters at z — 0, and did not examine their 
evolution. Considering the success of the enhanced star formation activity model, it is important 
to combine the ICM heating model with the cosmic star formation history. In doing that, we 
will have tighter constraints on the parameter space of e SN and e RG , and may have a link to 
the physical reasonable scenario beyond a general but parameterized modeling like our current 
approach. A future sample of clusters at high z selected by the Sunyaev-Zel'dovich effect 
(Sunyaev, Zel'dovich 1972) may provide another complementary piece of information on the 
thermal evolution of the ICM. We hope to report results on these important issues elsewhere 
in due course. 

This research was supported in part by the Grant-in-Aid for Scientific Research of 
JSPS (12640231, 14102004, 14740133, 15740157). Numerical computation was performed using 
computer facilities at the University of Tokyo supported by the Special Coordination Fund 
for Promoting Science and Technology, Ministry of Education, Culture, Sport, Science and 
Technology. 
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Appendix. Tracing thermal evolution of ICM in merging trees of dark matter 
halos 

In the present paper, we trace merger trees from z = to z = 30. In this case, however, 
the number of progenitors in one merging tree becomes progressively larger at higher redshifts, 
and it is not practically feasible. Thus we decide to prepare two different sets of merger trees 
at low redshifts (z = to z = 7) and at high redshifts (z = 7 to z = 30). 

We first construct N ens (=10) independent sets of merger tree realizations for each value 
of 12 different masses of M root starting from z min = to z max = 7 in N stcp = 750 timesteps. We 
also construct separately 20 sets of merger tree realizations starting from z min = 7 to z max = 30 
in iVstep = 500 timesteps. We use the number of timesteps in logarithmically equal redshift 
interval, A^ tep = 750 for the low z sets and N stcp = 500 for the high z sets. The masses of the 
root halo M root are equal to iq 1s +°- 2 ( 1 ~ - 5 ) M© (i = 1, . . . , 12) at z = for the low z sets, and 
2.0 x 10 7 " 1 "- 7 M Q (j = 1, . . . , 6) at z = 7 for the high z sets. The minimal mass of progenitors 
resolved in the merger tree M res is equal to M(T vir = 10 4 K), which corresponds to a mass 
of halos whose virial temperature is 10 4 K at each redshift. The parameters concerning the 
merging tree of halos is summarized in table 3. 



Table 3. Parameters for merger trees of dark halos. 



symbol 


adopted value 


physical meaning 




low redshift high redshift 




Mr 00 t 


1Q 13+0.2(i-0.5) Mq 2.0 X 10 7+J Mq 


mass of halo at z = z m j n 




(i = l,...,12) (i = l,...,6) 




Af ros 


Af(T vir = 10 4 K) 


minimal mass of progenitors resolved in each merger tree 


Astep 


750 500 


number of redshift bins (logarithmically equal interval) 


N cns 


10 20 


number of realizations of the merger tree 


^min 


7 


minimum redshift of the merger tree 


^max 


7 30 


maximum redshift of the merger tree 



Assuming that baryons at z = 30 have their gas temperature equal to the virial temper- 
ature of their individual host halo, both the halo merging history and the thermal evolution 
of gas are followed up to z = 7. At that redshift, each halo is assigned to another realization 
of the halo in the separate redshift merger trees using the interpolation of the mass in finding 
its counterpart. Then all the properties of the halo at high redshift tree is now transferred to 
its counterpart at the low redshift trees. Then we continue to follow its evolution until z = 
tracing the new merging tree. 

We checked the validity of the above interpolation method explicitly by constructing 
the merger trees from z = to z = 10, and tracing the evolution using each merger tree. The 
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comparison with the results based on statistically connecting two separate merger trees from 
z = to z = 7 and from z = 7 to z = 10 indicates that the properties of clusters at z = are 
almost indistinguishable. 
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